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Abstract 

The excitation energies of parahydrogen clusters have been systematically calculated by the 
diffusion Monte Carlo technique in steps of one molecule from 3 to 40 molecules. These clusters 
possess a very rich spectra, with angular momentum excitations arriving up to L = 13 for the 
heavier ones. No regular pattern can be guessed in terms of the angular momenta and the size of 
the cluster. Clusters with iV = 13 and 36 are characterized by a peak in the chemical potential and 
a large energy gap of the first excited level, which indicate the magical character of these clusters. 
From the calculated excitation energies the partition function has been obtained, thus allowing for 
an estimate of thermal effects. An enhanced production is predicted for cluster sizes N = 13, 31 
and 36, in agreement with experiment. 

PACS numbers: 67.40.Db, 36.40.-c, 61.46.Bc 



I. INTRODUCTION 



Small (pH.)„ clusters of parahydrogen have been produeed m a eryogenie free jet expan- 

sion and studied by Raman spectroscopy [1]. The Q{0) Raman line of the H2 monomer is 
shifted as the number of molecules in the cluster changes, thus providing a method to 
identify the cluster mass. The first seven resolved peaks next the monomer line have been 
assigned to clusters with = 2, . . . , 8 molecules. Although in that experiment the resolu- 
tion was not enough to resolve larger sizes, broad maxima were observed at A^ ~ 13 and 33, 
and perhaps 55, which have been interpreted as a propensity for geometric shell structures. 



Indeed, classical static 



BflQ 



J] and molecular dynamics 



results, based on a generic 



Lennard- Jones interaction potential, indicate that the expected structures for such clusters 
are the so-called Mackay icosahedra [6|, exhibiting some magic sizes (13, 33, 55...) related to 
the packing of molecules in closed icosahedral arrangements. A very complete discussion of 
classicalgeometrical patterns and their relation with the interaction features can be found 
in Ref. 7|] and references herein. However, quantum effects play a major role in pH2 clusters, 
like in their analogous helium droplets. The path integral Monte Carlo (PIMC) simulations 
of Sindzingre et al. [8| have shown indeed that the superfluid fraction in pH2 clusters with 13 
and 18 molecules become large at temperatures below T ~ 2 K. This prediction prompted 



both experimental 



91 and theoretical 



11 



12 



13| research of small clusters consisting of 



pH2 molecules surrounding an OCS chromophore, confirming the existence of a superfluid 
response. 

Quantum Monte Carlo (QMC) methods have been widely used in recent years as a the- 
oretical tool to study pH2 clusters. Several QMC techniques have been employed to cal- 
culate their pro per ties, namely variational Monte Carlo (VMC) and diffusion Monte Carlo 



(DMC) Q,n,[l6, 17, 18], path integral Monte Carlo (PIMC) [sl, Q, 20, 21, 3 



Monte Carlo (RMC) 



reptation 



23| and path integral ground state (PIGS) Monte Carlo [2J]. Only 



some specific values of the number of constituents A^ have been considered in the past, re- 
lated to the expected MacKay icosahedra structures. More recently, systematic calculations 
as a function of the number of molecules in the cluster have been performed. The ground 
state energies and the one-body densities of pHat clusters have been calculated by the DMC 
technique in steps of one molecule from N=3 to 5 [l7 |. by the PIGS technique from N=2 
to 20 y and by the PIMC from N=5 to 40 |2l|, I22J. The calculations show that (pH2)Ar 



2 



clusters exhibit a clear geometrical order, with the molecules occupying concentric spherical 
shells, which could be related to some polyhedric arrangement. The apparent incompati- 
bility between the large superfiuid fractions and the structured radial distribution densities 



has been recently clarified 



21 



22| by PIMC calculations which show that superfiuidity is 



localized at the surface of the clusters. 

Whereas up to ~ 22 these calculations are substantially in agreement, for heavier 
clusters there are noticeable differences between DMC and PIMC results, particularly for 

> 26. PIMC chemical potentials show very prominent peaks at N=26, 29, 34 and 39, in 

contrast with a smoother behavior obtained with DMC. These differences subsist even after 

I — I 

an improved DMC calculation ^], and seem to be related to thermal effects. According to 



Mezzacapo and Boninsegni 



21| they should be associated to a coexistence of solid-like and 



liquid-like phases, with a dominance of the latter at low T, as a result of both the zero-point 
motion and quantum permutation exchanges. However, thermal effects could manifest in 
enhanced stability thresholds at finite temperature, similarly to what has been observed in 
^He droplets 251]. Such droplets are definitely liquid-like and nevertheless they exhibit some 
magic numbers, which are not related to enhanced ground state binding energies at specific 
values of A^. They are instead stability thresholds, related to the cluster sizes at which 
excited levels cross the chemical potential curve and become stabilized [26]. An analysis of 



the possible existence of simi 
experimental results of Ref. 



ar enhanced stabilities in pH2 could be helpful to interpret the 
ll 

The purpose of the present article is to provide an accurate microscopic description of 
the ground state and the low-lying excited levels of (pH2)Ar clusters. In the past, the excited 
states of the N = 7 cluster were analyzed by McMahon et al. 27|]. We present here a 
systematic calculation of excitations with angular momentum from L = to 13 for clusters 
with A^ = 3, . . . , 40. For a given L-state the Schrodinger equation has been solved by means 
of stochastic DMC techniques using the pairwise potential determined by Buck et al. 28| . 
hereafter referred to as BHKOS. Another popular choice is the potential due to Silvera and 



Goldman 



29|, hereafter referred to as SC. Both potentials combine ah initio calculations 



with properties of the gas (or solid) as well as experimental information from collisions. 
The main difference among them is that SG contains a repulsive long-range term (cg/r^) to 



approximate the effective potential in a solid. As showed in Ref. 



the BHOKS potential 



provides more binding than the SG one, and the BHOKS chemical potential is slightly higher 



3 



than the SG one. 

The paper is organized as follows. In Sec. II, technical details of the DMC calculations are 
given for the ground state and low-lying excited states. The energy levels of (PH2) n clusters 
are presented and discussed in Sec. III. In Sec. IV the partition functions are obtained from 
the excitation energies, and they are applied to the calculation of the abundances in free jet 
expansions of parahydrogen clusters. A summary and some general conclusions are finally 
given in Sec. V. 



II. THE DIFFUSION MONTE CARLO METHOD 



The diffusion Monte Carlo (DMC) method 30| solves the imaginary-time Schrodinger 
equation for a function f{TZ,t) = $var('^)^('^, ^) which is the product of the importance- 
sampling wave function $var and the true ground-state wave function \E'. The variable TZ 
represents the set (ri,...,rjv) of 3A^ coordinates of the atoms which form the drop. 
The solution for / is obtained by constructing an approximate small-time Green's function 
G{71, 71', r) which serves to advance / by a small time step r, 

f{n',t + T) = j dnGin',n,T)f{n,t). (i) 

An initial set of walkers {7^i,7^2, • • - T^n^} is created to represent the function f(7l,0) at 
the initial time. After many repeated applications of the short-time approximate Green's 
function, all components of the starting wave function orthogonal to the ground state wave 
function disappear and the contributions of the other eigenstates go exponentially to zero as 
the time goes to infinity. The remaining set of walkers then provides a valid representation 
of the converged solution. In order to have a good statistical accuracy it is convenient to 
continue the application of the Green's function along a very large number of time steps and 
accumulate the resulting set of walkers. Then the mixed estimators of physically interesting 
quantities, such as the total energy, can be computed. Of course since the walkers are 
strongly correlated a block average is carried out in order to estimate the variance. 

Our ground state calculations are based on a Jastrow-like importance sampling wave 
function depending on just two parameters 

" ( 

2 \ b J p 



where Vij = |r.j — r^l is the distance between the pair Once the variational parameters 

were determined for each value of A^, DMC calculations were carried out with very small real- 
time steps (5 X 10^^ K^^ to 1 X 10^^ K^^, for the light and heavier clusters, respectively) and 
very long number of evaluation steps. Notice that the Bose symmetry of the above guiding 
wave function is not modified in the DMC process. The inclusion of triplet correlations 
in the importance sampling function leads to a noticeable improvement of the variational 



energies, as showed in Ref. 18j. However, the DMC energies are essentially the same as 
those obtained with only two-body correlations, with slightly more binding in the heavier 
clusters. The main difference lies in the reduction of the standard deviation by typically a 
factor of two. 

Importance sampling wave functions with a non-zero angular momentum may equally 
well be considered, and correspondingly the initial state will be an admixture of eigenstates 
of a specified angular momentum L. A convenient trial function with angular momentum L 
is the following 

N 

^i'^J{n)=J2^^{\^^-^cXyLL{0O^)} >< 'fvar(7^) , (3) 

i=l 

where Ui stands for the spherical angles of the vector r.^ — Rem which defines the position of 
the i-th particle with respect to the center-of-mass coordinate Rem = J2i ^i/^^ Ylm is a 
spherical harmonic. Note that the sum over all the constituent particles and the subtraction 
of the center of mass coordinate Rem are required to maintain the the Bose symmetry 
and the translational invariance, respectively. Again, the application of the imaginary-time 
Green's function will select the lowest eigenstate in the subspace of angular momentum 
L, the remaining contributions being exponentially suppressed as the time increases. The 
pref actor entering Eq. ([3]) is identical to zero when L = 1. In that case we have used 

N 

= I]Re {|r, - Rem|'rn(c^.)} x $,,,(7^) . (4) 

i=l 

A similar procedure to deal with L-states has been used in Ref. 
and in Ref. 26| for ^He droplets. 

Certainly, the trial function for L ^ states is not positive definite and consequently 
cannot be interpreted as aprobability distribution function. We have used the so-called 



27| for the iV = 7 cluster. 



fixed-node approximation 3l|] assuming that the nodal surfaces are fixed by the importance 



sampling function. In this approximation an exact energy for the lowest angular momen- 



turn L state is, in fact, not obtained. The improved variational wave function and the 
corresponding mixed estimator of the energy, which is derived, does, however, provide a 
variational upper bound. 

Besides the excited states with angular momentum L, one should also consider vibra- 
tional excitations characterized by a radial quantum number n. In short, let us represent 
each state by {n,L). With the previous trial wave function, the DMC method selects the 
lowest eigenstate (0, L) in the subspace of angular momentum L. An useful estimate of the 
energies of the average vibrational excitations {n > 1,L) in each L-subspace is obtained 



by using sum rules 32|, |33l], which only require the knowledge of the lowest L-state wave 



functions. Consider the exact lowest state for a given angular momentum L and the full 
set of eigenstates of this subspace in order of increasing energy, represented by \1/q^^ and 
En^^}, respectively, where n = 0,1 . . . is the vibrational quantum number. Let Q{Tl) 
be an arbitrary operator, assumed to be scalar under rotations, for which the sum rule of 
order p is given by 

where the sum extends over all eigenstates of the Hamiltonian with the exception of the 
lowest energy state (0, L) of angular momentum L. An upper bound to the energy of the 
first excited state of the subspace L is obtained as 



- E^'^ < ^7X7^ > (6) 



which involves the energy weighted and the non-energy weighted sum rules. Using the 
completeness relation of the set of functions {\E'^,^''}, they can be written as 



2 



= ^(^o'^lEl^^^n^o'V (8) 

j=l 

The computation of these expressions only requires a knowledge of the lowest state wave 
function of the angular momentum L subspace. 

We present here results for the first vibrational state of the subspace L = 0, employing a 
general monopolar operator 

p=l i<j 



The parameters Cp have been optimized so as to obtain the lowest upper bound to the first 
excited state. Notice that the excitations obtained are approximate because mixed matrix 
elements are used to calculate the required sum rules, so the strict variational character is 
lost. 



III. THE ENERGY LEVELS OF PARAHYDROGEN CLUSTERS 

Once calculated for each cluster the ground state E^^^ {N) and the excited state energies 
Eif'\N) and e[^\n), the excitation energies are defined as 

AElMN) = E^,'^\n) - e!,'\n) , (10) 
AEo{N) = Ef'\N) - (AT) . (11) 

They are plotted as a function of N in figures [l]l3l For the sake of clarity, we have displayed 
them separately for each value of L in figures [Tll2] and the full spectra is showed in figure [31 
The chemical potential or dissociation energy 

MiV)=Er(iV-l)-E(°)(iV), (12) 

defined as the difference between the ground state energies of neighbor clusters, has also 
been plotted in all cases. Indeed the chemical potential defines the energetic region in which 
excited states are stable. Only stab^ 

As previously discussed in Refs. 
ical potential is the enhanced stability of the cluster = 13, which confirms its magical 
character. A mild peak appears also at = 36, but one cannot exclude that it could simply 
be a statistical fluctuation. An additional confirmation of its magicity will be given later 



e states have been plotted in figures [Tll3l 

ItI . [isl , the most prominent result regarding the chem- 



on. In contrast, PIMC calculations of Refs. 



2l| and 



22l |. at T=l and 0.5 K respectively, 



show an enhanced stability at = 13, 26, 29, 34, and 39. 

The main characteristics of the spectra displayed in figures [1] and [2] can be summarized 
as follows. All clusters exhibit stable excitations for L = 2, 3 and at A^ = 3 and onwards. 
The next excited level with L = 4 starts to be bound at A^ = 4, the following levels with 
L = 5 and 1 appear at A^ = 6. The next L = 6, . . . levels appear at regularly increasing 
size thresholds. In the considered size range the highest stable excited level corresponds to 
L = 13, which appears at A^ = 31 and onwards. The quadrupolar L = 2 excitation is the 



lowest one for all sizes, except at = 26, 28, 29 and 37, for which the octupolar L = 3 
state lies at a slightly lower energy. The excitations with L > 7 do not cross each other: the 
higher the angular momentum the higher the excitation energy. The vibrational excitation 
energies of the state {n = 1,L = 0) present a smooth behavior as a function of N, apart from 
a mild minimum at = 13. We have not pushed our calculations to the vibrational states 
(1, L 7^ 0). Presumably, they would also present a monotonic behavior with respect to their 
corresponding (0, L) states. In contrast, the excitation energies of the states {n = 0, L ^ 0) 
have a non-monotonic behavior with A^. 

The = 13 cluster is very peculiar. Its excitation energies present a prominent peak 
for the states with L = 1 — 5, a prominent minimum for L = 6, and a pronounced change 
of trend for L = 7, 8 and 9. This cluster has no stable excited states for L > 9. No other 
so enhanced bumps appear up to A^ = 40. Such a behavior is related to the truly magic 
character of this size. Its enhanced stability results in a repulsion of the excitations, similarly 
to what has been observed in classical Lennard- Jones calculations % IS^]- The prominent 
minimum found for L = 6 at A^ = 13 means that this value of the angular momentum is also 
special. A similar phenomenon was observed in the study of doped parahydrogen clusters 
C0@(pll2) N by reptation QMC calculations 23|], where the time-correlation functions of the 
various multipole moments in the rotating-axes frame have been computed. A persistent 
crystal structure for A^ = 12 and a persistent time correlation of the L = 6 multipole has 
been observed by these authors. Their explanation for these facts was based in the property 
that the only non-null multipolar moments for the A^ = 13 cluster with icosahedral symmetry 
are L = 6 and L = 10. 

Our results are compared in Table [T] with those of McMahon et al. 27| for N = 7 
which, to the best of our knowledge, are the only available calculations of excited states of 
parahydrogen clusters. These authors have performed a DMC calculation using the BHKOS 
interaction and a trial wave function for L-excited states analogous to our Eq. ([3]). It can 
be seen in the table that both calculations are in substantial agreement within at most two 
standard deviations. 

In figure [3] are collected our calculated DMC excitation energies, omitting error bars for 
the sake of clarity. The obtained spectrum for each (pH2)Ar cluster is complex, in the sense 
that it does not resemble any of the simple models of excitations. The states with L = 2 
are the lowest lying, apart from the few cases mentioned before. It is worth mentioning 



8 



the presence of relevant energy gaps at = 13 and = 36, in coincidence with peaks in 
the chemical potential. Both facts indicate a magical character of these clusters. Note that 

r 

analogous energy gaps have been determined in classical Lennard- Jones clusters [4|. 

There is almost no information, neither experimental nor theoretical, of the spectrum 
of clusters weakly bound by van der Waals like interactions. An indirect determination 

3^]. These interesting 



35 



36 



of the (pH2)2 dimer excitations has been reported in Refs. 
experiments refer to the absorption of infrared radiation by parahydrogen gas at 20 K, 
slightly above the triple point. Just focusing in the region around the Qi{0) line there has 
been observed the absence of absorption at exactly the frequency of the free Qi{0) line, and 
the presence of two satellites at energies above and below the free parahydrogen frequency, 
interpreted as a signal of the dimer excited states. In practice, absorption experiments 
provide an indirect way of detecting the two-body ground state pH2— pH2 and measuring 
its spectrum. 

In the case of ^He droplets, an indirect determination of their excitation energies results 
from the temperature-dependent analysis of the production abundances of such clusters 
in ultrasonic expansion of pressurized gas 26|]. Sudden jumps have been observed in the 
quantum partition function at cluster sizes which can just accommodate one more additional 
stable excitation. The nice agreement between the experimental magic sizes and theoretical 
stability thresholds provides a confirmation of the calculated energy levels of these clusters. 

In principle one could expect some parallelism between the spectra of ^He droplets and 
those made of parahydrogen: both constituents are bosons and the interaction is rather 
similar. However, our calculations rule out such an analogy, as helium clusters present 
regularly spaced spectra. In contrast, no simple pattern can be guessed for the spectra of 
parahydrogen clusters. Only the excitation energy of the L = level has a smooth behavior 
with the number of pll2 molecules, with a small inflection at = 13. Its value is roughly 
constant near 20 K, quite far away from the dissociation limit, in contrast with the analogous 
level in ^He droplets which is close to the dissociation limit. 



The (n = 0, L > 1) of ^He excitations correspond to sur 



like oscillations with radial frequencies given by 



38 



39, 



40| 



ace modes, described by ripplon 



.1/2 



These frequencies were corrected by Tamura 



L(L-l)(L + 2) 



(13) 



4l| to include finite size effects and the 



corrected expression reproduces very well the calculated DMC spectra of helium clus- 
ters j26|. Certainly, the smooth behavior related to the angular momentum dependence 
[L{L - 1)(L + l)]i/2 does not emerge in our calculated spectra of pH2 clusters. 

The excitation energy of the translationally invariant forbidden 1 mode, which requires 
a nodal surface in the radial coordinate space, has a quite high energy compared with the 
allowed levels (2,3,...), like in the case of helium droplets. Nevertheless, its energy is again 
far away from the dissociation limit, in contrast with the case of ^He clusters. With respect 
to the variation with the number of molecules, it is noticeably smoother than the levels of 
low angular momentum (up to L = 4) but reveals a prominent bump near = 13. 

At the light of our calculations we infer a coexistence of liquid and crystal structure 
near = 13, the liquid part being the responsable of the L = 2,3, ... excitations, and the 
crystal part the responsable of the L = 6 excitation. The region around A^ = 13 where 
the energy of the L = 6 mode is decreasing up to A^ = 13 and afterwards increasing up to 
a plateau will correspond to partial icosahedral symmetry. We have analyzed the known 
crystal structures related to Lennard- Jones clusters reported in the Cambridge Cluster Data 
Base 42[, by computing the multipolar moments. There are three significative cases worth 
to be mentioned. Firstly, the A^ = 19 cluster, for which the multipolar moments L = 3 and 
L = 4 are null and for which we obtain a very low L = 2 excitation and much higher L = 3 
(a peak) and L = 4. Secondly, the A^ = 26 case, whose first non-null multipole is L = 3 and 
appears in our calculations with a minimum for this excitation. Finally, the A^ = 38 cluster, 
whose first non-null multipole is L = 4 appearing again in our calculations a valley around 
this value of the number of constituents. In conclusion, even if we do not find a simple 
or familiar model to describe the excitations of parahydrogen clusters, some aspects of the 
spectrum are qualitatively described by assuming a coexistence of liquid and solid phases. 



IV. THERMAL EFFECTS AND MAGIC NUMBERS 

The knowledge of the full spectra allows us to analyze different thermal effects by means 
of the partition function. In particular, we consider in this section the expectation value of 
the temperature-dependent energy E{T) and the effect of excited states on the production 
abundances of clusters in free jet expansion of pressurized gas. 
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The partition function of each cluster is expressed as the product 



N 



(14) 



of the center-of-mass (Z^) and internal {C,n) partition functions. The former can be written 

as 

^^=Kr^ • (15) 



The internal partition function is given in terms of bound state energy levels of the state 
j = (n, L) of the N-molecule cluster by 



in = ^gj{N) exp 



kT 



(16) 



where the degeneracy factor for the state with the radial quantum number n and angular 
momentum L is given by (/j = 2L + 1. 

By using the notation of section UTTl Eqs. fllOHllj) . the internal partition function may be 
written in a computationally convenient form 

where AEl{N) are the excitation energies for angular momentum L, and the statistical 
factor has been written in the explicit form (2L + 1). 



N 



exp 



l + ^(2L + l)exp 



kT 



(17) 



The internal energy of a given cluster at a temperature T is given by 



1 + El(2^ + l)AEL(Ar)e-^^^(^)/'=^ 



(18) 



1 + El(2^ + l)e-^^^W/^^ 
Sizable thermal effects on the energy values will start to appear at temperatures close to 
the energy of the first excited state. Thus, at T = 1 the energy changes are minimal, with a 
maximum of 0.4 K for N ^ 20, related to the minimum of AE2 (see fig. [1]). The energy of 
the other clusters remains practically unaltered, and consequently the chemical potential at 
T = 1 has practically the same shape as at T = 0. This observation is interesting in relation 
with PIMC calculations, which have been carried out at T = 0.5 K 22! and T = 1 K 21|. 
The series of peaks observed in the PIMC chemical potential as a function of the number 
of molecules is not reproduced by the present extension of DMC to T = 1. Although 
different interactions have been used in these calculations, namely BHKOS in DMC and SG 
in PIMC, we do not expect this as the origin of the differences. Indeed, DMC results using 
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both interactions [18| give substantially the same chemical potential shape, apart from the 
mentioned mild peak at = 36, not visible with the SG interaction. 

On the other hand, the effect of the temperature on the mechanism by which clusters 
are formed in the free jet expansion turns out to be very interesting. The analysis which 



follows is similar to the one carried out for Helium droplets |25l . |26| which explained the ex- 
perimentally observed anomalously high production rates of clusters with a specific number 
of constituents, the anomalies being related to the excitation spectrum of the clusters and 
not to a magical character. 

Let us assume that cluster growth is dominated by the chemical equilibrium reaction 

ipH2)N-i +pH2 + X^ {pH2)n + X , (19) 

where X is a spectator particle required to fulfill the energy and momentum conservation 
laws. The related equilibrium constant is independent of the spectator particle and is given 
by a quotient of partition functions 

Km = . (20) 

The equilibrium constant is dominated by the quotient Stv/Sat-i, which is conveniently 
written as 



Neither the chemical potential nor the excitation energies are smooth functions of A^, as 
has been show in figure [31 One can expect that the above quotient will reflect such a non 
monotonic behavior as jumps at particular values of A^. The quotient of partition functions 
Sjv/Sat.i is displayed in figure [Has a function of A^ for several values of temperature, lower 
than the source temperatures of Ref. about 40 K. As these ratios span several orders of 
magnitude as T varies, we have arbitrarily given the value 1 to the case A^ = 40. 

Clearly pronounced peaks are observed at A^ = 13, 31 and 36, and less pronounced at 
A^ = 26. Other peaks exist, like A^ = 19 and A^ = 29, but they only appear at specific 
values of T. The pronounced peaks give an explanation of the experimentally determined 
production peaks at A^ = 13 and 33. It should be mentioned that the relation between 
the Raman shift and the number of molecules of the clusters in the work of Ref. [l| is 
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an extrapolation from small values of N, and the position of the peak at = 33 may be 
imprecise. It is worth noticing that the last value is half way between the theoretical peaks 
at iV = 31 and N = 36 related to an enhanced production. 



V. CONCLUSIONS 

In this work we have determined the excited-state energy levels of parahydrogen clusters 
containing up to 40 molecules by means of a DMC calculation based on an importance- 
sampling function with a nonzero angular momentum L plus the fixed-node approximation. 
A variational upper bound has been obtained for the lowest-lying energy of a given L-state. 
Bound excited states have been found for states (n = 0, L = 1 — 13) and {n = 1, L = 0), 
the size threshold for binding increasing as L increases. It turns out that the spectra are 
complex, in the sense that no simple pattern can be guessed to fit them as a function of 
the angular momentum of the state and the size number. In particular, and contrary to 
what happens in the analogous system of helium droplets, the liquid drop model formula is 
excluded. 

The = 13 and A^ = 36 clusters are magic, as indicated by the existence of peaks in 
the chemical potential in coincidence with energy gaps for their first low-lying excited state, 
relative to their neighbors. 

The spectra of weakly bound van der Waals clusters have not yet been experimentally 



measured. Only an indirect determination 25|, l26(] has resulted from the analysis of enhanced 
production of clusters of ^He in terms of the successive stabilization of excited levels at 
growing values of the number of constituents of the clusters. The analysis of the partition 
function of parahydrogen clusters has revealed the existence of production enhancements at 
A^ = 13, 31 and 36, in nice agreement with the experiment l|]. The effect is similar to the 
case of helium clusters, but the nature is very different. Instead of the stability threshold 
reason of helium clusters, here the enhanced production is basically related to the existence of 
energy gaps in the excitations at the mentioned values of the number of constituents. There 
may be some corrections due to not including the {n = 1,L ^ 0) levels in the partition 
function. We expect, however, these effects to be small, as presumably such vibrational 
excitations would be some 20 K above the non- vibrational (n = 0, L) corresponding level. 
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f. Q 



TABLE I: Cluster with N = 7. DMC results of Ref. [27l | are compared with present results. All 
energies are in K 

McMahon et al. Present results 

-98.89(15) -98.76(1) 

H 25.90(15) 25.37(2) 

AEa 8.3(3) 7.92(13) 

/^E^ 13.2(2) 12.92(17) 

A^5 15.4(2) 14.74(18) 
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FIG. 1: (Color online) Excitation energies (in K) of (pH2)iv clusters as a function of the number N 
of constituents for states with L = 1 — 6. The chemical potential n{N) is also displayed to indicate 
the stability limit. 
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FIG. 2: (Color online) Excitation energies (in K) of (pH2)Ar clusters as a function of the number 
N of constituents for states with L = 7 — 13 and the vibrational state with L = 0. The chemical 
potential n{N) (green dots) is also displayed to indicate the stability limit. 
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FIG. 3: (Color online) Excitation energies (in K) of (pH2)Ar clusters as a function of the number of 
constituents. The chemical potential n{N) (green dots) is also displayed to indicate the stability 
limit. Levels with L = 2 and are plotted in red and blue, respectively. 
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FIG. 4: (Color online) The ratio of the partition functions Hjv/Stv-i is plotted versus N for several 
values of temperature: 1, 5, 10, 15, and 20 K. The scale has been arbitrarily fixed at = 40. 
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